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More than others, arctic ecosystems are affected by consequences of global climate changes. The herbivorous plays 
numerous roles both in Scandinavian natural and cultural landscapes (Forbes et ah, 2007). Wild reindeer (Rangifer 
tarandus L.) herds in Hardangervidda plateau (Norway) constitute one of the isolated populations along Fennoscandia 
mountain range. The study aims to understand temporal and spatial variability of intra- and inter-annual home ranges 
extent and geophysical properties. We then characterize phenological variability with Corine Land Cover ecological 
habitat assessment and bi-monthly NDVI index (MODIS 13Q1, 250 m.). Thirdly, we test relationships between rein- 
deer’s estimated densities and geophysical factors. All along the study, a Python toolbox (“GRiD”) has been mounted 
and refined to fit with biogeographical expectancies. The toolbox let user’s choice of inputs and facilitate then the 
gathering of raster datasets with given spatial extent of clipping and resolution. The grid generation and cells extraction 
gives one tabular output, allowing then to easily compute complex geostatistical analysis with regular spreadsheets. 
Results are based on reindeer’s home ranges, associated extent (MODIS tile) and spatial resolution (250m). Spatial 
mismatch of 0.6 % has been found between ecological habitat when comparing raw (100m 2 ) and new dataset (250m 2 ). 
Inter-annual home ranges analysis describes differences between inter-seasonal migrations (early spring, end of the 
summer) and calving or capitalizing times. For intra-annual home ranges, significant correlations have been found 
between reindeer’s estimated densities and both altitudes and phenology. GRiD performance and biogeographical 
results suggests 1) to enhance geometric accuracy 2) better examine links between estimated densities and NDVI. 


Keywords: Raster datasets compiling, biogeography, Rangifer tarandus L., home ranges, kernel densities estimation, 
ecological habitat, phenology 


1. Introduction 

In the information age, « big data » constitutes a large part of the future of ecology and biogeography (Hampton et 
al., 2013). Transdisciplinary expertise for socio-ecological issues are nowadays relying on multi -factorial and multi- 
scale analysis for environmental surveys. Massive datasets are easily computed and processed, thanks to tremendous 
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computer hardware advances and trends for horizontal interactions between experts, citizens and institutions (Kam- 
batla et ah, 2014). That said, the era of data-intensive science needs some specific tools to collect and centralize the 
numerous and heterogeneous datasets related to a studied phenomenon, such as wildlife and herding issues, interacting 
with cultural and natural systems. Multi-temporal scalar information is another issue, long-term strong trends (climate 
change, biodiversity loss) interacting with inter-annual variability. Within the scope of environmental impact assess- 
ment, and particularly in the field of animal ecology, an equal spatial resolution of biotic and abiotic factors is required 
to describe, explain and forecast how a population is supposed to move in a limited range of space, according to future 
bioclimate and biocenoces interactions. 

Arctic and subarctic ecosystems are particularly affected by the consequences of global climate change (Uboni et 
al., 2016). Rangifer tarandus tarandus L. constitutes a key-component species, and its role in both natural and cultural 
sphere is increasingly important. With a wide circumpolar distribution, the capital-breeder and long-distance migra- 
tory ungulate moves through taiga and tundra biomes. Its particular geographical distribution and its ecological habi- 
tats selection are relevant to ecosystem structure and functioning. The ability of Rangifer to cope with its food needs 
and spatial requirements may represent a good opportunity to increase arctic and subarctic socio-ecosystems resilience 
to rapid effects of climate change (Post and Pedersen, 2008). At a global scale, populations of wild reindeer are clas- 
sified as vulnerable by the International Union for Conservation of Nature (UICN), and currently red-listed by the 
institution. Warmer winter temperatures, rapid changes in the water cycle, increasing landscape fragmentation and 
human disturbances contribute to deplete populations, biocenoces communities and individuals (Weladji et al., 2006). 

On the Hardangervidda plateau (Norway), the subpopulation of Rangifer tarandus tarandus L. is one of the last 
remaining wild reindeer population in Scandinavia. The population stock has been reduced sevenfold between early 
1970s and today (Uboni et al., 2016). With Rondane and Snohetta sub -populations, those three biogeographical isles 
could represent a genetic resource contributing to Scandinavian Rangifer resiliency to climate change and related 
cultural and natural landscapes. As such, it is important to improve our understanding of the capital -breeder ungulate 
migrations during the vegetation growth period and the critical times in its life cycle: response to spring snow melt, 
calving period and “fat accumulation” period in late summer (Klein, 1990; Skarin, 2008). Our hypothesis are the 
following: Rangifer tarandus tarandus L. has a seasonal behavior, choosing different habitats according to the key 
moment of its biological cycle during the early spring and the growing season (Skarin, 2008); strong estimated densi- 
ties of Rangifer tarandus tarandus L. depend on the bioclimatic conditions during its migration, and on geographical 
features such as topography and habitats (Klein, 1990). 

Management of merged datasets, here by spatializing relationships between biocenoces and biotopes, may improve 
GIS toolbox as well as our understanding of Rangifer tarandus interactions. The study aims to present a GRiD-toolbox, 
combining some pre-existent tools in current versions of ArcGIS and QGIS software, applied to a multi-scalar reindeer 
home ranges study. 


2. Material and methods 


2.1 Material 

Several free-access datasets, fully described in Table 1, are compiled, merged and analyzed in the study. Down- 
loading datasets is the first step (url access in Table 1). The first basis is the reindeer GPS -tracks database used by 
Cagnacci et al. (2015) to assess and characterize most European’s largest herbivorous seasonal residence times, asso- 
ciated home ranges and identifying migrational structures. 7 wild female reindeers have been equipped and surveyed 
using GPS-radio collar. The second dataset is ASTER-GDEM, providing altitudes and other topographical infor- 
mations (slopes, expositions). The third dataset involves average annual temperatures, interpolated from the Global 
Historical Climatological Network at high resolution (Hijmans et al., 2005). Corine Land Cover map is the fourth 
dataset. The Corine program assess ecological habitats and land uses at European scale, with remote sensing and 
photo-interpretation observations. The downloaded datasets are available for 1990, 2000, 2006 and 2012. In our case, 
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the 2006 edition is used in order to suit temporally with the locations of wild reindeer. The last dataset is represented 
by MODIS 13Q1 time series. Bi-monthly temporal acquisition has been chosen, from the beginning of 2000 to the 
end of 2015. Dates of imagery acquisition is expressed in Julian day. The acquired tiles extent to entire Scandinavia, 
with a spatial resolution of 250 meters. Since the native dataset is in*.HDR format, Normalized Difference Vegetation 
Index (NDVI) band has been chosen, and converted into GeoTiff format. 

Different softwares are used in the study, depending on processing steps. First of all, we used GIS software, such 
as ArcMap 10.2 and Qgis 2.18 for both preprocessing and processing. Then, R-studio and Ade-Habitat-HR (Home 
Ranges, Callenge, 2006) package allow us to calculate reindeer home -ranges, using Kernel Utilization Density tech- 
nique, on a bi-monthly basis. Then home ranges have been exported in separate GIS files (*.SHP). These home ranges 
and associated quantification of reindeer’s densities constitute the key- input for the generation of new datasets. Reg- 
ular spreadsheets, such as Open Calc 4. 1 .3 or Office Excel 2016 versions are used all along the study. Both R (package 
Rcmdr) and Xlstat can be used for complex statistical calculation. For this specific study, Xlstat has been privileged. 

Table 1 . Studied datasets and associated metadata 


Dataset 

type 

Reindeer GPS -tracks 

DEM 

Annual tem- 
perature 

Ecological habitats 

Vegetation phenol- 
ogy 

Sampled 

variable 

“KUDrefvalue” 

“alti”, “expo”, 
“slopes” 

“tmoy” 

Corine Land 

Cover habitats 

MODIS 13Q1- 
NDVI, 2007-2009 
growing season 

Format 

*.txt 

*.tif 

*bil 

*.tif 

*.hdr 

Data 

source 

http : //datadryad. org/re- 
source/doi: 10.506 l/dryad.rg0v3 

http://www.eea.eu- 

ropa.eu/ 

www.world- 

clim.org 

http://www.eea.eu- 

ropa.eu/ 

https ://earthex- 
plorer.usgs.gov/ 

Period 

available 

3/20/2007-2/15/2010 

- 

1960-1990 

2006 

2000-2015 

Geo- 

graphical 

extent 

Hardangervidda plateau (Norway); 
up to 100 km 2 

World, tiles 

World 

Europe 

World, tiles 

Temporal 

resolution 

3 hours 

- 

Annual 

1990-2000-2006- 

2012 

Bi-monthly 

Spatial 

resolution 

-meters 

30 meters 

900 meters 

100 meters 

250 meters 

Unit of 

measure- 
ment/type 
of varia- 
ble 

- 

Meters, degrees 

Celsius 

Categories 

NDVI index 

Geodesic 

system 

WGS84-UTM 32N 

ETRS 89 

WGS 84 

ETRS 89 

Polar sinusoidal 


2.2 Methods: defining spatial units , extracting values and statistical analysis 

The general approach requires to define the spatial and temporal unit to compare and to analyze, in this case bi- 
monthly home ranges during the growing season (April to August). Grid extent will be generated at the regional scale 
and related cells used for the extraction will depend on bi-monthly home ranges distribution. The spatial and temporal 
resolution is forced here on the MODIS- 13Q1 accuracy, but this can be specified in GRiD-toolbox parameters. 


2.2.1 Monitoring wild reindeer ( Rangifer tarandus L.) home ranges in Hardangervidda Plateau (Norway) 
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At the scale of the entire population equipped with GPS-collars (seven individuals), the GPS-tracks has been sorted 
out by bimonthly period and exported in *.CSV format. Data have been sorted in two packages. The first sorting 
aggregates bi-monthly GPS-tracks during the whole available period (2007-2010) and represent the seasonal variabil- 
ity of the growing period (two bi-monthly home ranges by month, from April to August, totaling ten files). The second 
sorting focuses on the inter-annual variability of the ten bi-monthly home ranges during three years (2007, 2008 and 
2009), giving 30 files. 

The GPS-tracks are then ready to estimate the reindeers’ utilization distribution. Using the R-package ADE Habi- 
tatHR (Callenge, 2006), we compute the Kernel Utilization Distribution (KUD) for each fortnight. The bivariate nor- 
mal kernel technique is specified, as well as the smoothing parameter calculated with « ad-hoc » method (Worton, 
1989; Skarin et al., 2008). Grid size (-100 m 2 ) and extent (h-value for estimated densities; 1km 2 ) are also specified to 
fit with reindeer’s large scale migration properties. The resulting convex hull has been vectorized, with 95% of the 
kernel home-range extent as boundary. Then datasets are exported, respectively into GeoTiff raster grid for Kernel 
Utilization Distribution densities, and into shapefiles for the 95% convex hull. Finally, a grid is generated, intersected 
following reindeer KUD’ shapes and geophysical descriptors systematically sampled. 


2.2.2 GRiD - toolbox: ArcMap modelbuilding and Python programming (Qgis) 


The figure 1 summarizes the processing chain used for GRiD ( Grid Raster Information Dataset) tool-box creation. 
Two prototypes of GRiD are available, even if improvements are still on progress. Here the ArcMap Modelbuilder is 
preferred, but Qgis GRiD toolbox is ready to use for the third and last step shown in figure 1 . According to the goals 
of the study, the first three steps are flexible in the process chain. 

A key step is to define a common projection in order to minimize geometrical error between datasets. The WGS84- 
UTM 32N projection has been chosen, initially used for locations of wild reindeers. After projecting every dataset, 
we generated a regional grid, based on MODIS 13Q1 tile extent and spatial resolution (250m). Each cell is theoreti- 
cally supposed to suit with MODIS 13Q1 coarse grid, and thus values of NDVI. The step 3 converts features (cells) 
to points located at the core of each cell. The resulting mesh of points still has a regional extent, such as the basic grid 
obtained from step 2. The generated bi-monthly home ranges of wild reindeers are put into the chain, in order to 
intersect the regional mesh of points. As output, each intersect processing produces one shapefile and one spatial unit 
for home range analysis. 


Definir une 
projection 
Defining a 


* 


common 

projection 


• From each raster 
dataset projection... 

• ...to home ranges 
one's 



Grid generation 


1 At regional scale 
(MODIS tile extent 
and spatial 
resolution) 

> Each cell will be later 
filled with 
geographical 


4 



core of each grid's cell 


1 point = one 250 
meter MODIS 13Q1 
value (ideally) 


1 


4 



! JIL 4 

am 

Intersect 

• Crossing bi-monthly 
home ranges (annual 
or pluri-annual) 

• And the general mesh 
of points 


tin <<_1 )*!?*>*• 


m jm» ann 1W v 

m nts« rows mjjsi k 

m 296J67 » 6»7 7/MS7 7. 

599 29656* 29656* (.164*6 71 


Sample 

• For every geophysical 
parameter 

• 1 CSV file = 1 
reindeer's home 
range (annual or 
pluri-annual) 


Fig. 1. General process line for GRiD toolbox, using ArcGIS Model-builder (credits: R. Courault, 2017) 

The final step is crucial: raster values are sampled and resulting attribute tables exported as a *.CSV file. Each 
generated CSV matrix blend spatial unit of analysis (here bi-monthly home ranges of reindeers) with different geo- 
physical datasets we are willing to statistically sum up or establish relationships (columns). Each observation (line) is 
a pixel contained into a specific pre-calculated home range. These geo statistical observations are now characterized 
with one specific value both for biotic and abiotic factors, and are ready to analyze. 


2.2.3 Statistical processing: grid validation and redundant information removal for explanatory factors 
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2 . 2 . 3. 1 Checking correspondences between raw and sampled datasets 

We quantify the possibility to lose spatial information, by upscaling geophysical datasets at MODIS - 13Q1 spatial 
resolution (-250m). In the study, the comparison will be made between raw Corine Land Cover dataset (e.g. 100m of 
spatial resolution) and new Corine dataset extracted using GRiD toolbox (250m). The analysis involves the calculation 
of particular and total error (%) between raw and new dataset. A threshold of 5% of total error has been retained. 


2. 2. 3. 2 Removing geophysical factors which are spatially autocorrelated 

The goal is here about deleting geophysical factors which could be redundant in the analysis between reindeer’s 
home ranges and reindeer’s kernel densities within bi-monthly home ranges. At the scale of one of the most extent 
home range (2 nd half of April), continuous values (annual temperatures, altitudes, exposition, slopes, NDVI) are 
crossed using Pearson’s correlation. A threshold of 0.8 is retained for the significance of a relationship, allowing then 
to remove redundant factors. 


2.2.4 Statistical processing at intra-annual scale of analysis 


2.2.4. 1 Characterizing the seasonal variability of the chosen ecological habitats 

For the characterization of the inter- seasonal variability of the choice of ecological habitats, we summed every 
pixel belonging to each habitat class within the Kernel Utilization Distribution area. Adding up pixels allow us to 
compare habitat’s areas (sum of pixel’s absolute values, or percentage calculation) between seasonal home ranges, 
aggregated during the whole period. The statistical summaries and bar charts for every home range are then summed 
up by selecting the first and the second most frequent ecological habitat by bi-monthly reindeer’s home range. A Khi 2 - 
test is then computed to verify the independence between the fortnight and the category of land cover (number of 
pixels). 


2. 2.4.2 Inter-annual variability of bioclimatic conditions: choice of particularly different years for spring migration, 
calving area and “fat accumulation” area 

When extracting MODIS-NDVI time series for annual home ranges, care was taken to sample the three Julian date 
acquisition relevant to bi-monthly home range. It gives three NDVI time series (columns) by bi-monthly home range. 
For example, when extracting geophysical raster based on 1 st half of April, the 97 th Julian day of 2007, 2008 and 2009 
have been added to the 1 st half of April home range. We then check the independency of these three inter-annual 
statistical distribution, using Kruskal -Wallis independency tests. Based upon the median of rank generated distribu- 
tions, two different years in regards of NDVI inter-annual distributions have been selected. 


2.2.5 Statistical processing at inter-annual scale: specific bi-monthly home ranges 


2.2.5. 1 Testing relationships between biotic and abiotic factors 

The selected contextual bi-monthly home ranges correspond to three major moment in reindeer’s annual biology 
(descent to summer pastures = 1 st half of April, calving time = 2 nd half of May, fat accumulation time = 2 nd half of 
July). For those particular home ranges, we tested the relationships between contextual estimated densities and other 
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extracted factors (altitudes, exposition, slopes or proxies of primary chlorophyll production such as NDVI). Spearman 
correlation test is preferred in this case. According to Spearman table of significance, with a confidence threshold of 
95% and a number of observations >100, the given threshold for significance is above 0.197. 


2. 2. 5. 2 Differences in habitat choices 

We compare the NDVI values distribution and the reindeer’s densities probabilities distributions belonging to the 
same home range. Our aim is to verify the similarities or dissimilarities between these two kinds of distributions. 
Kruskal- Wallis independence tests has been conducted using habitat categories as sub-samples. Ecological habitat 
categories divide both densities estimation and NDVI values as many sub -samples as ecological classes the observa- 
tions belong to. 


3. Results and suggestions 


3.1 GRiD toolbox: sufficient accuracy 


Table. 2. Accumulated error, in percentages, between raw Corine Land Cover dataset (spatial resolution =100 m 2 ) and extracted Corine 
Land Cover using GRiD toolbox (spatial resolution = 250 m 2 ) for the 2 nd half of April 


% of total area for each present class of habitat for raw Corine Land Cover dataset 


Dis- 

contin- 

uous 

urban 

fabric 

(2) 

Sport 

and lei- 
sure fa- 
cilities 

(11) 

Non-ir- 

rigated 

arable 

land 

(12) 

Complex 

cultiva- 
tion pat- 
terns 

(20) 

Land 

princi- 
pally oc- 
cupied by 

agricul- 
ture (21) 

Broad- 

leaved 

forest 

(23) 

Conif- 

erous 

forest 

(24) 

Mixed 

forest 

(25) 

Moors 

and 

heath- 

land 

(27) 

Transi- 

tional 

wood- 

land 

shrub 

(29) 

Beaches, 

dunes, 

sands (31) 

Sparsely 

vegetated 

areas (32) 

Peat 

bogs 

(36) 

Water 

courses 

(40) 

Water 

bod- 

ies 

(41) 

To- 

tal 

0.04 

0.05 

0.07 

0.30 

0.77 

16.51 

9.49 

5.80 

38.96 

0.55 

0.24 

15.14 

5.45 

0.02 

6.60 

100 

% of total area for each present class of habitat, extracted Corine Land Cover dataset, using GRiD 

0.04 

0.04 

0.08 

0.32 

0.75 

16.40 

9.41 

5.84 

39.16 

0.57 

0.24 

15.14 

5.44 

0.04 

6.54 

100 

Accumulated error for each present class of habitat and total of overall accumulated error 

0 

0.01 

0.01 

0.01 

0.02 

0.11 

0.08 

0.03 

0.2 

0.03 

0 

0 

0 

0.02 

0.07 

0.6 


Table 2 compares the percentages of Corine Land Cover ecological habitats (1 st line, table 2) according to the 
spatial resolution, in one of the most extent home range (2 nd half of April, e.g. 17236 cells; -4309 km 2 ). Similar 
percentages can be found within coarse dataset (3 rd line) compared to the percentages for GRiD -extracted dataset. The 
last line subtracts the difference between GRiD -extracted percentages to coarse Corine Land Cover. Each accounted 
error has been added up to the total, giving a total error about 0.59 %. Close to 99.4% of representativeness, the GRiD- 
toolbox extraction shows good results, in particular for Corine Land Cover ecological habitats. 

That said, some improvements could reinforce spatial accuracy of GRiD -toolbox. Firstly, the quantification of error 
between specific spatial resolution of coarse datasets and inputted spatial resolution using GRiD-toolbox has to be 
more thoroughly assessed. It could be interesting to repeat the same test for other datasets, such as raw ASTER-GDEM 
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(spatial resolution e.g. 30m) and extracted values (250m). At the scale of the analysis (here the 2 nd most extent rein- 
deer’s home range), it is important to notice that the number of observations (e.g. pixels) can vary according to the 
considered spatial unit. Since the number of geostatistical observations is very different between home ranges, check- 
ing spatial unit representativeness for each dataset would allow us to have a good overview of total error for a given 
spatial resolution as GRiD toolbox input parameter. 

The second improvement is related to geometrical bias, probably generated by the chosen map projection (UTM- 
32N). Indeed, Mercator projection is more suited with equatorial region, and tends to elongate surfaces when dealing 
with polar or subpolar regions. Working with an official map projection clipped by zone (EUREF89-NTM in our case) 
will likely improve overall results for GRiD -toolbox and wild reindeer distribution analysis. Another parameter 
(choosing cartographic projection) could be added to the toolbox and automatically process the projection homogeni- 
zation between the numerous datasets we are willing to extract. 


3.2 Intra-annual scale: towards a selection of biotic and abiotic parameters 


3.2.1 Auto-correlated and redundant data 

The auto-correlated factors are annual mean temperatures values (1950-1990) and altitudes. This can be explained 
by the altitudinal gradient, in particular between the high altitudes and low annual temperatures in the north of the 
studied area, low altitudes and higher temperatures in the south (Hardangervidda plateau). At the scale of the second 
most extent home range, the correlation is strong and significant (r = +0.83; p<0.0001). Then, it is statistically possible 
to express altitudes distribution as a good proxy for annual temperature distribution. That said, it could be interesting 
to reprocess the calculation with a more recent dataset (1980-2010 for example) and more accurate temperatures var- 
iables (minimal and maximal temperatures). 


3.2.2 Ecological habitats and intra-annual (e.g. inter-seasonal) variability of home ranges 

In absolute values of pixels, intermediary seasonal home ranges (early spring, late summer) have the greatest over- 
all area, compared to calving period and full-season home ranges. For example, the most extent home range (1 st half 
of April) counts 19736 geostatistical observations (-4480 km 2 ) when the smallest one (2 nd half of June) is about 2262 
pixels (-565 km 2 ). Indirectly, this fact is shown in the figure 2. Percentages of the two most frequent ecological 
habitats by fortnight and home ranges during the whole period are here displayed. Colors represent the type of eco- 
logical habitat, with a total number of 6. For both April and August fortnights, a relatively balanced distribution is 
noticed compared to June. Moors and heathland areas are widely represented, with no less than 35% as minima for 
every home range total areas. Water bodies are represented for the 1 st July home range (-10%), since any internal 
differentiation has been made within bi-monthly home ranges. 

The Khi 2 -test (performed on the number of pixels of the 9 more represented habitats for the ten fortnights) shows 
a strong relationship between these variables (observed value: 9873.29, critical value: 92.80, p<0.0001). Some habitats 
are over-represented during early spring, particularly sparsely vegetated area. During spring and summer, moors and 
heathland and mixed forests are over-represented. During late summer, coniferous forests and transitional woodlands 
are over-represented. 

Thus, due to the wider extent of reindeer’s home ranges during intermediary seasons, habitat diversity appears to 
be higher when females start to go towards the calving areas. The same can be said during the summer, with mixed 
forests and later coniferous habitats emerging after the phenological peak of broad-leaved forests followed by the one 
of moors. 
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% of area for the two most frequent habitat, by bi-monthly home range 



1st April 2nd April 1st May 2nd May 1st June 2nd June 1st July 2nd July 1st August 2nd August 


■ Moors and heathlands Sparsely vegetated area ■ Broad-leaved forests 

■ Mixed forests ■ Water bodies ■ Coniferous forests 

Fig. 2. Comparison between reindeer’s fortnight home ranges along vegetational season (April to August) and the first two major frequent 
CLC habitat classes by home range (sources: Corine Land Cover ecological habitat, European Environmental Agency; AdeHabitat pack- 
age, Callenge, 2006; GRiD-toolbox; credits: R. Courault, 2017) 


3.2.3 NDVI dissimilarities and selecting specific home ranges 

The Kruskal- Wallis test is used to verify if statistical distributions are significantly different. This is the case for 
the three annual home ranges tested (spring migration: 1 st half of April; calving area: 2 nd half of May; fat accumulation 
area: 2 nd half of July), totaling 9 distributions. The NDVI values (2007-2008 and 2009) are statistically different 
(p<0.0001). Between the three years, it is now possible to determine what are the « cold » and « warm » years accord- 
ing to NDVI onset, distribution and central parameter. These chosen years will be studied for the inter-annual home 
ranges. The specific home range describing migration toward summer pastures (1 st half of April) under “cold” condi- 
tions is 2008, with +0.1 1 as average of the NDVI. 2009 is the “warm” year with +0.32 as average of NDVI. The home 
range describing calving period (2d fortnight of May) under “cold” conditions is 2008 (average: +0.45) and 2009 
“warm” (+0.56). The specific home range period describing the fat accumulation period under “cold” conditions is 
2009 (+0.66), whereas 2007 is “warm” with +0.75. 


3.3 Inter-annual scale of home ranges: what parameters are related to reindeer estimated 
densities? 


3.3.1 Spring as a robust marker for estimated densities of wild reindeers 

Table 3 reports Spearman’s rho correlations between estimated reindeer’s densities (Kernel Utilization Distribution 
technique) for selected years and relevant explicative factors. Correlations exceeding a threshold of 0.2 (p-value 
<0.0001) are shown here. High altitudes are likely associated with high values of reindeer densities during cold years 
1 st fortnight of April and 2 nd of May, year 2008). High-altitude plant communities, such as bryophytes and terricolous 
lichens could be preferred, reindeers awaiting seasonal meltdown and vascular plants onset (Klein, 1990). Conversely, 
“warm” conditions are most likely to correlate with estimated densities of reindeers and NDVI. In that case, correla- 
tions describe both negative (Spearman’ s-p =-0.268, p-value < 0.0001) for April 2009 and positive Spearman’s cor- 
relation (+0.262, p-value < 0.0001). The negative correlation (1 st fortnight of April 2008) depicts low densities related 
to high NDVI values, whereas the positive correlation (2 nd fortnight of May 2009) depicts high reindeer densities with 
high NDVI values. Explanations can be found when looking at opportunistic reindeer feeding strategies, in particular 
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in spring when annual herbaceous just start growing (Skarin et ah, 2008). No clear correlation for “fat accumulation” 
home ranges (2 nd fortnight of July) has been found. Here, the relative abundance of biomass, shortly after the pheno- 
logical peak, can explain the lack of relationship. 

Table 3. Spearman’s rank correlation matrix between densities estimation by annual home ranges and explicative parameters. Only Spear- 
man’s rho > 0.2 and statistical p-value <0.0001 are shown here 


Estimation of utilization distribution by 
annual home ranges/ Explicative factors 

Altitudes 

NDVI 

1st fortnight of April 2008 (cold) 

0.374 

NA 

1st fortnight of April 2009 (warm) 

NA 

-0.268 

2nd fortnight of May 2008 (cold) 

0.230 

NA 

2nd fortnight of May 2009 (warm) 

NA 

0.262 


3.3.2 Strong inter-annual variability in estimated densities and phenology between the inter-annual home 
ranges 


In order to explore more thor- 
oughly the relationship between 
reindeer estimated densities and 
NDVI values, Kruskal- Wallis tests 
have been computed. Statistical dis- 
tributions have been clustered by 
sub-samples, e.g. ecological habitat 
classes. Every test calculated for the 
9 specific home ranges has rejected 
the dependence hypothesis between 
the geostatistical distributions of es- 
timated reindeer densities and 
NDVI values (9 tests for the 9 spe- 
cific home ranges). Such variability 
is displayed in figure 3. The set of 
maps shows likeness in pattern dis- 
tribution of estimated densities 
shapes (ethological/individual be- 
havior, calving area), but also dis- 
similarities (environmental factors 
depending on contextual year, 
spring migration). 



Fig. 3. Set of maps showing reindeer estimated densities using Kernel Utilization Distribution technique, in a color scale from blue to red. 
a. 1 st fortnight of April 2008; b. 1 st fortnight of April 2009; c. 2 nd fortnight of May 2008; d. 2 nd fortnight of May 2009; e. 2 nd fortnight of 
July 2007; f. 2 nd fortnight of July 2009 (Base-map sources: ESRI Dataglobe, Credits: R. Courault, 2017) 
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4. Conclusions 

Mounting GRiD toolbox has been incrementally encouraged by the needs of automation processing for the survey of 
wild reindeer migration and densities in Hardangervidda Plateau (Norway). The wide range spatial behavior of rein- 
deer is an interesting case study to apply this method, more adapted to large territories than simple grid multi -data 
analysis previously realized by other authors in smaller areas (Jolivet 2014, Bortolamiol et ah, 2016). Indeed, by 
choosing the core point of each grid cell as the common reference to cross all the datasets, the process is realized 
within a reasonable computation time while avoiding important bias according to the preliminary test realized in our 
study. Improvements still have to be made to ensure GRiD geometrical reliability, cartographic projection homogeni- 
zation and minimizing error when extracting large and various raster datasets. 

Applying GRiD toolbox for wild reindeer survey during the three growth periods simplified pre -treatments by 
automatizing grid and creation of mesh of points; and finally datasets fusion into one attribute table. Based on geosta- 
tistical observations (pixels) it allowed us to compute several complex statistical analyses, not usually present in GIS 
software. Based on Open Access datasets, the study firstly shows interesting links between reindeer intra-annual dis- 
tribution and ecological habitat diversity; secondly relationships between inter-annual variability of estimated densi- 
ties and abiotic factors, validating thus our second and third hypothesis. 
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